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Level and wavefunction statistics have been studied for two dimensional clusters of the square 
lattice in the presence of random magnetic fluxes. Fluxes traversing lattice plaquettes are distributed 
uniformly between — i<E>o and i$o with $0 the flux quantum. All considered statistics start close to 
the corresponding Wigner-Dyson distribution for small system sizes and monotonically move towards 
Poisson statistics as the cluster size increases. Scaling is quite rapid for states close to the band edges 
but really difficult to observe for states well within the band. Localization properties are discussed 
considering two different scenarios. Experimental measurement of one of the considered statistics 
-wavefunction statistics seems the most promising one- could discern between both possibilities. 
A real version of the previous model, i.e., a system that is invariant under time reversal, has been 
studied concurrently to get coincidences and differences with the Hermitian model. 



I. INTRODUCTION 

There is a general belief based on Anderson localiza- 
tion theory Q that all states of two-dimensional (2D) 
systems are localized in the absence of a magnetic field. 
The situation is not so clear when time-reversal symme- 
try is destroyed by a magnetic field. The motion of a sin- 
gle particle in a 2D random magnetic field is attracting 
strong interest since it is related both to the half-filled 
Quantum Hall Effect [||,[| and the slave-boson descrip- 
tion of high T c superconductors The existence of 
recent experiments measuring transport properties in a 
static random magnetic field adds considerable interest 
to this subject ||. 

While perturbative renormalization group calculations 
show that all states are localized |||7| , the eventual pres- 
ence of an extra term could give rise to a Kosterlitz- 
Thouless transition from the localized phase to a phase 
with power law correlations Numerical simulations 
have used a variety of forms: the iterative strip method 
has been used by Sugiyama and Nagaosa Avishai et 
al. Jig ], Liu et al. jllj], and Kalmeyer et al. anal- 
ysis of participation ratio has been employed by Pryor 
and Zee 12 and Kalmeyer and Zhang tails in the 
density of states have been discussed by Gavazzi et al. 
|l3f and Barelli et al. JTJ], a network model has been 
introduced by Lee and Chalker jl5) , the diffusion of elec- 
trons has been studied by solving the time-dependent 
Schrodinger equation [fl6|| , and finally, Hall conductivity 
has been used as a way for determining energy regions 
where extended states dominated in the thermodynamic 
limit fLTjj . Loosely speaking, states look extended in fi- 
nite samples as soon as the edge of the band is left and 
a definite answer is quite difficult. Therefore, the most 
extended conclusion is that a mobility edge separates lo- 
calized states from extended or critical states. On the 
contrary, Sugiyama and Nagaosa |)| conclude that all 
states are localized in random magnetic fields based on 



the non-existence of a second branch in the scaling func- 
tion whereas the same conclusion is reached by Lee and 
Chalker jl5| based on finite values of the localization 
length obtained in the semiclassical limit described by 
their model. 

In this work, I present numerical evidence showing 
that, in the presence of random fluxes, states scale 
towards random (Poisson) statistics everywhere in the 
band. Several standard statistics (wavefunctions ampli- 
tude statistics, nearest neighbor spacing statistics and 
number variance) are used in conjunction with other cri- 
teria (for example, spatial extension of wavefunctions as 
given by the participation ratio) to show a systematic 
tendency towards Poisson statistics (the one that char- 
acterizes systems with a spectrum formed by exponen- 
tially localized states) as the size of the random sam- 
ples increases. Scaling is quite rapid near the band edge 
where states "look" exponentially localized whereas it 
seems to follow a logarithmic law near the band center 
where states "look" fractal jL8|. Nevertheless, since my 
analysis is based on statistical magnitudes, the eventual 
existence of a countable set of extended states can not be 
disproved. In other words, the results presented in the 
next Sections prove that states scale towards localization 
on average. 

Recent results regarding level statistics can be found in 
the literature jl^jist] . Most of them deal with the char- 
acterization of a new universal distribution function that 
would be appropriate at the metal-insulator transition. 
As the results presented in the rest of the paper show, 
signs suggesting the existence of a critical energy (mo- 
bility edge) at which scaling changes sign have not been 
found. Scaling flows in the same direction in all energy 
regions. The paper is organized as follows. Section II de- 
fines the model that has been numerically solved in the 
paper, Section III gives the main results grouped under 
three subjects (Density of states, wavefunctions and level 
statistics) and finally, main conclusions are collected in 
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the last Section. 

II. LATTICE MODEL 

The Hamiltonian describing random gauge fields on a 
L x L cluster of the square lattice is : 

H=-tJ2^"'c^i' > (!) 

<w> 

where cj creates an electron on site I, I and I' are nearest- 
neighbor sites, and —t is the hopping energy (hereafter 
t = 1 is chosen). The flux through a given loop S on 
the lattice is <f>s = X)<z i>>es measured in units of 
flux quantum <J> = ^f- Link fields satisfy ipw = — 4>i'i- 
Although some checks have been performed for the Meiss- 
ner phase |l3[ (link fields uniformly, randomly distributed 
in the interval [—5,5]), the bulk of the calculations have 
been done for the Debye phase in which uncorrelated 
fluxes are randomly selected from the interval [—hi §]• 
Notice that in any case, the situation of maximum pos- 
sible disorder is considered except when the density of 
states is analysed as a function of disorder. Model (|l|) 
will be called RMF (Random Magnetic Fluxes or Ran- 
dom Magnetic Field) model in the rest of the paper. 

In order to recover standard results, I have always run 
checks for a real version of Hamiltonian (|l]) . To this end, 
link field values are restricted to and 1/2 and chosen 
randomly. The model describes a square lattice in which 
the hopping integral is constant in absolute value but 
takes random signs. This model will be shortly men- 
tioned as RHS (Random Hopping Signs) model. Stan- 
dard localization theory should apply to this model. Both 
models are equivalent for a particular choice of random 
fluxes; when <j>ui is always equal to 4 but has random 
sign, hopping integrals take random ±2 complex values. 
Thanks to the bipartite character of the square lattice, 
this Hamiltonian can be made real: it suffices a canoni- 
cal transformation of all operators on one of the sublat- 
tices from fi; to i&i. In this way, this particular complex 
Hamiltonian is described by the same matrices as the 
RHS model. This case provides an overlap between real 
random Hamiltonians that according to standard theory 
show complete localization and a Hamiltonian with com- 
plex eigenfunctions. 

III. RESULTS 
A. Density of states 

One of the easiest ways allowable to characterize the 
effect of disorder is the density of states (DOS): 

N(E) = ^S(E-E a ) (2) 



where the sum runs over all eigenstates of the system. 
For arbitrary disorder {4>w} the eigenvalues of Hamilto- 
nian ([l]) lie in the interval [—4, 4]. Moreover, there is also 
an exact symmetry on a bipartite lattice like the square 
lattice. For every eigenfunction of energy E there is an 
eigenfunction of energy —E whose amplitudes have op- 
posite sign on one of the sublattices. As a consequence, 
the density of states is symmetric about E = 0. The 
same symmetry forces the existence of an eigenstate of 
zero energy in clusters of an odd number of sites. 
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FIG. 1. Evolution of the density of states as the inter- 
val [— 0max, 0max] from which random values of the flux 
are selected increases. From top to bottom: ideal system, 
0max = g, \, § and |. The DOS of a Bethe lattice of coordi- 
nation equal to four has been plotted (continuous line) in the 
bottom panel to allow comparison with the RMF result. 

Generally, an average over disorder realizations is done 
in order to get a smooth result. Averaging is particularly 
important when Anderson model is studied for large dis- 
order because band width increases and the number of 
states per energy unit diminish. Since the disorder given 
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by Hamiltonian ([j]) systematically decreases the width of 
the band, the study of only one large sample suffices to 
get nice results. Fig.(|f|) shows the evolution of the den- 
sity of states obtained for a 200 x 200 system as disorder 
(random fluxes distributed in the interval [— max , </>max]) 
increases. A small amount of disorder destroys the log- 
arithmic singularity at the center of the regular system 
and shrinks somewhat the band. Peaks appear near the 
band edges due to some kind of accumulation of states 
due to the band skrinkage. Finally, a smooth curve is 
obtained for the maximum possible amount of disorder: 
random fluxes between — \ and \. In this case, the den- 
sity of states of a Bethe lattice of coordination equal to 
four: 



N(E) 



2 V12 - E 2 
7T 16 -E 2 



(3) 



gives a very good approximation to the random result (it 
is difficult to distinguish between histogram and continu- 
ous curve in the bottom panel of Fig.(|l|)). Although the 
characterization and deep understanding of the peak ap- 
pearing close to the edges at small disorder deserves fur- 
ther study, the rest of the paper is devoted exclusively to 
the case of maximum disorder (bottom panel of Fig. ([!])). 



matrix is not too large (a typical workstation can deal 
with a 200 x 200 cluster, i.e., with a 40000 x 40000 matrix 
if only eigenvalues are calculated but only with a 64 x 64 
system (a matrix of dimension equal to 4096) if eigen- 
vectors should be obtained). On the other hand, inverse 
iteration |2lJ allows the calculation of particular eigen- 
vectors. The algorithm starts with some random values 
for the components of the wavefunction ?/>o of energy E 
and proceeds as follows: 



= {E-H)- 1 ^ 
V> 2 = {E- H)- 1 ^! 



(4) 



Convergence is reached when |^ — \ is smaller than 
some tolerance e (e = 10 -8 works fine in double preci- 
sion). Energy should be updated after some number of 
cycles (about 10): 

Ei =< ipi\H\ipi > . 



FIG. 2. Typical "extended" wavefunction of a 250 x 250 
cluster and energy E w — 7r/50. 



B. Wavefunctions 

Wavefunctions can be obtained using standard algo- 
rithms for Hermitian band matrices [|20|| if the size of the 



FIG. 3. Typical "localized" wavefunction of a 250 x 250 
cluster and energy E w — 7r/50. 

Taking advantage of the band structure of the Hamil- 
tonian, eigenstates for clusters as large as 300 x 300 can 
be obtained using standard linear algebra subroutines. 
Sometimes larger systems are accesible: for example, 
since E = is always an eigenenergy of Hamiltonian 
([j]) when the number of sites is odd, total number of it- 
erations is greatly reduced and eigenvectors for 399 x 399 
clusters can be obtained. The computational thresholds 
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that have been given refer to complex Hermitian matri- 
ces. Somewhat larger systems can be studied for the real 
version of model ([!]) . 

Wavefunctions produced by Hamiltonian (|l|) present 
a great richness. We get bona fide exponentially local- 
ized wavefunctions for energies close to the band edges, 
almost localized wavefunctions (i.e., wavefunctions that 
show non vanishing amplitudes in a small part of the 
cluster) in any part of the energy spectrum, quasi one- 
dimensional eigenstates in many cases and also, quasi 
extended wavefuntions in the main part of the spec- 
trum. Figs.(H) and(||) show typical results at cluster sizes 
at which only inverse iteration allows the calculation of 
specific eigenstates Wavefunctions are represented 

through a 256 values grey scale in which darker means 
larger values of \ipi\ 2 being ipi the wavefunction ampli- 
tude on site I. Although a large collection of such wave- 
functions have been compiled, it is clear that the simple 
display of them does not help to the analysis of local- 
ization. Therefore, I have use standard tools to analyse 
their extension and the way they behave as the size of 
the sytem is increased. 



1. Participation ratio 

A good measure of the spatial extension of an eigen- 
state is given by the participation ratio: 



LxL 



1=1 



(5) 



where ipi is the amplitude of wavefunction on site I. P 
can be understood as the number of sites covered by the 
wavefunction and reaches a maximum finite value for lo- 
calized states. In order to avoid systematic errors due 
to the use of a biased definition, I have also tried alter- 
native definitions of the spatial extension covered by a 
wavefunction. One is given in Casati and Molinari work 



(LxL \ 
-^|V/| 2 log(|V/| 2 ) 
1=1 J 



(6) 



a definition of the number of sites asociated to a wave- 
function that emphasizes small amplitudes and gives con- 
sequently larger values than the previous definition. Still 
another definition is based in the basic understanding of 
\ipi\ 2 as the probability of finding the particle at site I. 
We proceed in two steps: (i) weights are ordered in 
decreasing magnitude and (ii) , S is defined as the number 
of ordered weights that have to be added to reach some 
percentage of the total weight 1 of the state. I have tried 
different percentages to learn that this simplest defini- 
tion gives results for S that are practically identical to C 
when 99% of the total weight is considered. Using 90% 



as the cutoff, results closely follow those obtained using 
the standard P definition. This knowledge gave me fur- 
ther confidence in the reliability of the usual definitions 
of number of sites covered by a wavefunction. 
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FIG. 4. Number of sites for all states of seven samples of 
64 x 64 size. Top: standard definition of the participation 
ratio(|^). Bottom: Logarithmic definition 

Fig.(Q) and Fig. (j^) show both P and C for a collection 
of the largest clusters for which the complete set of eigen- 
ergies and eigenvectors can be obtained using my limited 
computing facilities. Schwarz algorithm for symmetric 
or Hermitian band matrices has been used jjjcj]. These 
results show that spatial extension is small near band 
edges but increases quickly towards the band center. In 
fact, the behavior is really sharp for the complex RMF 
model. One can be tempted to speak about some kind 
of pseudo-mobility edge at these energies (~ ±3) where 
spatial extent sharply increases. In fact these energies 
coincide with the values given in the literature as the 
limits of a reliable calculation of the localization length 
(see, for example, the paper by Sugiyama and Nagaosa 
||). Moreover, while states well within the band show 
a steady power law increase in their spatial extension, 
states close to the edges do not extend as the size of the 
system is increased. 
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FIG. 5. Number of sites for all states of two samples of 
58 x 58 size. Top: standard definition of the participation 
ratio(|^). Bottom: Logarithmic definition 

The scaling behavior of the spatial extension of states 
has been carefully analysed at particular energies. In this 
analysis inverse iteration is used and therefore, eigen- 
states for much larger systems are available. Fig.(^) 
collects results for RHS model at several energies while 
Fig.(^) shows spatial extension results for RMF model 
just at the band center. Let me begin the discussion with 
some comments about the results displayed in Fig.(||). 
Several energies covering different parts of the spectrum 
have been chosen. While it is clear that the growth of the 
spatial extension is bound in the lower part of the band 
(E = —3.4 and perhaps E = —3), results at energies well 
within the band show eigenstates that extend more and 
more as the system size increases. Of course, the validity 
of this conclusion is limited by the sizes that are reach- 
able by a numerical study (the largest cluster analysed 
in this case is 210 x 210). Nevertheless, a quadratic fit to 
the log-log plot shows a small negative curvature that in- 
dicates that after a region in which < P >~ L a (a < 2), 
a region of true exponential localization of wavefunctions 
could follow. Again a word of caution is in order at this 
point: the inference of the asymptotic behavior of cover- 
age from the analysis of a limited range of sizes is really 



a hard question in any numerical study. Therefore, the 
only conclusion that seems to be sure is that coverage in- 
creases following a power law with an exponent smaller 
than 2 at the most. 
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FIG. 6. Spatial extension of eigenstates produced by the 
real version of Hamiltonian (jl|) at different energies. Circles 
give the value of log(P) of the state that is closest to the 
selected energy for each sample of a set of realizations of dis- 
order. L is the cluster side. Solid lines are quadratic fits to 
the log-log values. 

A similar analysis is displayed in Fig. (]7|) for the Hermi- 
tian Hamiltonian. Participation ratios for states at the 
band center are given for increasing cluster sizes. In that 
case, the maximum size allowed by our computing capa- 
bilities is 399 x 399. Even for this system size, the number 
of sites covered by eigenstates at the band center is fol- 
lowing on average a law close to < P >~ L 16 . The cur- 
vature of a quadratic fit to the log-log data shows again 
negative -although really small- curvature (solid line of 
Fig. ([?])). So, we are forced to conclude that numerical 
methods based in the study of the scaling properties of 
individual states are not enough to show exponential lo- 
calization in 2D random gauge fields. 



2. Fractal dimension 



Fractal dimension has been defined as: 
rflog(P) 

Ct = 

dlog(L) 



(7) 
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It characterizes energies or system sizes for which the 
number of sites covered by eigenstates increases follow- 
ing a power law (Number of sites oc L°). It has been 
calculated for the whole spectrum using energy averaged 
values for the number of sites and two different system 
sizes to calculate the derivative. In other words, we have: 



a 



A log < P > 
A log < L > 



(8) 
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FIG. 7. Spatial extension of eigenstates produced by 
Hamiltonian Eq.(|l]) at the band center. Diamonds correspond 
to different realizations of disorder for several sizes (L x L 
cluster). Solid line is a quadratic fit to the log-log values. 

Fig.(||) gives our numerical results for both models. 
Sizes L = 32, 64 have been used for RHS model whereas 
sizes L = 24, 48 did the work for the more difficult RMF 
model. It can be seen that although the overall behav- 
ior is similar for both models, the increase to « 1.8 is 
so sharp for the RMF model that the idea of some criti- 
cal energy separating exponentially localized states from 
states following a power law decay could be considered. 

The fact that fractal dimension for averages of the cov- 
erage P is below 1.8, together with the more detailed 
results at particular energies shown in the previous sub- 
subsection, proves overall localization at least in a power 
law way: 

<P> 



This behavior is in striking contrast with the result that 
applies to wavefunctions of the GUE: 



< P > 
N 



where the total number of sites N plays the role of L 2 in 
the previous expression. For the same reason, localiza- 
tion is also present in the simplified RHS model for which 
fractal dimension is always less than 1.4. Actually, states 
are less extended for this model. 
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FIG. 8. Fractal dimension of both models (RHS: grey line; 
RMF: black line) as a function of the position of states within 
the band. 



3. Wavefunction amplitude statistics 

Wavefunction statistics measures the spatial varia- 
tion of eigenfunctions over the system under study. 
Crystalline wavefunction have constant amplitude 1/N 
whereas wavefunction of disordered or chaotic systems 
show large fluctuations. The wavefunction statistics of 
the canonical Gaussian ensembles is given by the corre- 
sponding Porter-Thomas distribution [f23| . In particular, 
we have 



/(*) = 



: exp(-t/2) 



(9) 



for the Gaussian Orthogonal Ensemble, where t is the 
squared wavefunction amplitude divided by its mean 
value, i.e., (t — Nip(r) 2 ) being TV the number of 
sites. Graphical representation is improved if variable 
x = log(t) is used instead of t. We get: 



L 2 



V27T I ^ 



(10) 



G 



On the other hand, we have: 
f(t) = exp(-i) 



(11) 



for the Gaussian Unitary Ensemble, where t is the 
squared wavefunction amplitude divided by its mean 
value, i.e., (t = Nip(r)*tp(r)). Doing the same variable 
transformation as before, we get: 



g(x) = exp(a; — exp(a;)) 
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model while Fig. (11) gives RMF results. In both cases, 
the same clear tendency of diverging from Porter-Thomas 
statistics as size increases, can be seen. Nevertheless, the 
distance to a typical statistics produced by exponentially 
localized wavefunctions is still large (thick solid line of 
Fig. (pi]) shows wavefunction statistics at the band edge). 
Notice that exponential localization implies a statistics 
that shifts towards left (smaller weights are more proba- 
ble) as the size of the system increases beyond the local- 
ization length. 
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FIG. 9. Wavefunction statistics of RMF model. Distribu- 
tion of wavefunction weight for i) an average over all states of 
256 8 x 8 samples (dotted line), ii) an average over all states 
of 64 16 x 16 samples (solid line), iii) an average over all 
states of sixteen 32 x 32 samples (grey thick solid line) and, 
iv) an average over all states of two 58 x 58 samples (black 
thick solid line). GUE Porter-Thomas distribution is given 
for comparison (dashed line). 

Deviations from these canonical distributions have 
been recently calculated by Falko and Efetov |24]. Their 
work applies to pre-localized states in disordered con- 
ductors. Numerically, it is a simple task getting devia- 
tions from universal behavior. Fig.(|^) shows the rapid 
deviation of RMF wavefunction statistics from the corre- 
sponding Porter- Thomas distribution (Eq.([l2"|)) as soon 
as eigenfunctions of medium-size matrices are considered. 
In this calculation, all states -regardingless their energy- 
have been used to compute the probability density. Close 
to the maximum capability of my computer (58 x 58 clus- 
ters), it can be seen that the ratio of evolution of wave- 
function statistics with size is decreasing. 

The divergence of wavefunction statistics from the cor- 
responding Porter-Thomas distributions shown in the 
last Figure is confirmed by the calculation of wavefunc- 
tions of much larger clusters close to a fixed energy (I 
have used inverse iteration starting at E = — 7r/50 try- 
ing to get unbiased eigenstates that could exist at spe- 
cial energies). Fig.(|To|) shows selected results for RHS 



FIG. 10. Wavefunction statistics of RHS model. Distri- 
bution of wavefunction weight for i) an average over all states 
of sixteen 16 x 16 samples (dotted line), ii) an average over 
all states of four 32 x 32 samples (solid line), iii) an aver- 
age over fifteen states of energy E ~ — 7r/50 obtained for 
different 200 x 200 clusters (grey thick solid line) and, iv) 
an average over ten states of energy E ~ — 7r/50 obtained 
for different 300 x 300 clusters (black thick solid line). GOE 
Porter-Thomas distribution is given for comparison (dashed 
line) . 

C. Level statistics 

1. Nearest neighbor spacing statistics 

Level statistics has been recently used to characterize 
the properties of spectra near the mobility edge p9| . It 
is well known that the distribution of nearest neighbor 
spacings is of the Wigner-Dyson type |25[] when parti- 
cles move through the whole slightly disordered sample 
and changes to Poisson when states become exponen- 
tially localized. Simplifying the argument, we can say 
that only states localized in different spatial regions do 
not interact through the Hamiltonian and are, therefore, 
allowed to lie at the same energy. According to this the- 
ory, level statistics should move towards Poisson as the 
size of the system increases whenever states are local- 
ized in the thermodynamic limit whereas it should move 
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towards the corresponding Wigner-Dyson distribution if 
states are extended in the infinite system [^6| . 
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FIG. 11. Wavefunction statistics of RMF model. Distri- 
bution of wavefunction weight for i) an average over all states 
of sixteen 16 x 16 samples (dotted line), ii) an average over 
forty states of energy E w — 7r/50 obtained for different 50 x 50 
clusters (solid line), iii) an average over six states of energy 
E « — 7r/50 obtained for different 300 x 300 clusters (grey 
thick solid line) and, iv) an average over twenty states of en- 
ergy E « —3.4 obtained for different 50 x 50 clusters (black 
thick solid line). GUE Porter-Thomas distribution is given 
for comparison (dashed line). 

Fig.(|l~2|) shows the results of a standard analysis of 
nearest neighbor spacings produced by the RMF model. 
The statistics of small clusters (144 sites) closely follows 
Wigner surmise whereas statistics of the largest system 
that we can analyse (40000 sites) shows a small but clear 
deviation from GUE statistics. On the other hand, larger 
scaling features can be observed if the level spacing dis- 
tribution obtained at the lower part of the spectrum of 
the RHS model is analysed (see Fig.(|l3|)). Level statistics 
for small sizes (empty circles) reasonably follows Wigner- 
Dyson distribution (dashed line) whereas larger sizes give 
rise to a distribution (filled circles) that is approaching 
Poisson distribution (dotted line). 

2. Scaling of the variance of the spacing distribution 

Since the visual inspection of nearest neighbor spacing 
distribution is not a quantitative way of analysing their 
evolution with system size, I have selected the variance 
of the distribution as the scaling magnitude. Matrices 
belonging to the GOE (GUE) give a distribution of vari- 
ance equal to 0.286 (0.180) P7II , whereas the variance of 
the Poisson distribution is 1. So, a continuous increase 
of the width of the level spacing distribution (measured 
here by its variance) should denote the path extended- 



to-localized followed by wavefunctions. Large samples of 
nearest neighbor spacings has been obtained by direct di- 
agonalization of Hamiltonian matrices obtained for ran- 
dom realizations of disorder. The total number of lev- 
els has been maintained almost constant (~ 200000 for 
RHS model and ~ 480000 for RMF model) as the size 
of the square cluster has been increased from 12 x 12 to 
200 x 200. 
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SPACING 

FIG. 12. Top panel: Distribution of nearest neighbor spac- 
ings obtained for 4096 randomly generated 12 x 12 clus- 
ters (histogram). Bottom panel: Distribution of nearest 
neighbor spacings obtained for eleven randomly generated 
200 x 200 clusters (histogram). Poisson (dashed line) and 
Wigner-Dyson (solid line) distributions are shown for com- 
parison purposes. 

Firstly, I present the results obtained for the variance 
of the whole sets of spacings. They are given in Fig.jT^) 
for both models together with the results obtained for 
rectangular clusters (these are the only results for this 
cluster shape presented in the paper). While scaling is 
obviously towards localization for 2 x L systems, asym- 
totic values for squares (RHS model) and circles (RMF 
model) cannot be inferred in an unique way. A good fit 
is obtained using the following model curve: 

<r = A + C x /L + C 2 /\og{L)+ 
+C 3 /(L log(L)) + CV L 2 + C 5 /log(L) 2 , 

where a means variance and L is the cluster size. Unfor- 
tunately, while results for RHS point towards a value 
A = 1 in a three parameter fit (ideal localization of 
all eigenstates) , RMF results for A are not so well de- 
fined (the number of parameters makes appreciable dif- 
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ferences) but point to an A value smaller than 1 although 
a five parameter fit with fixed A = 1 also goes through all 
circles. If this overall fit could be taken quite seriously, it 
would mean that some percentage of the total number of 
states (about fifty percent since A ss 0.6) does not follow 
a random statistics even in the thermodynamic limit. I 
will come back to the discussion of this point after the 
presentation of all results. 

1.0 
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FIG. 13. Distribution of nearest neighbor spacings in the 
lower part of the band produced by the real version of Hamil- 
tonian (^) for two different system sizes: open circles corre- 
spond to an ensemble of 12 x 12 clusters whereas filled circles 
show the distribution for an ensemble of 200 x 200 clusters. 
Poisson and Wigner GOE distributions are shown for com- 
parison. 

Secondly, an analysis that takes into account the en- 
ergy region from which level spacings originate has been 
done. Level spacings have been grouped in seven sets 
corresponding to seven energy regions and the variance 
of the corresponding distributions has been calculated. 
Fig. (|l5|) shows the flow towards Poisson statistics fol- 
lowed by the level statistics within the different energy 
regions. Scaling in the lower part of the band shows a 
rapid tendency to localization while states near the band 
center show a much slower ratio. This behavior easily ex- 
plains why localization is easily detected near the band 
edges but is almost invisible near the band center (see 
Pig-®)- 

Numerical analysis is still more difficult for the Her- 
mitian RMF model since scaling towards localization is 
much slower for this model. Remember that the number 
of level spacings is large and almost constant (~ 480000) 
for all considered system sizes. In this case, the variance 
should start close to the value of 0.180 that corresponds 
to the nearest neighbor level spacing distribution given by 



matrices of the Gaussian Unitary Ensemble [^7j . Fig. ( |l6| ) 
shows the evolution of the variance of nearest neighbor 
level spacings within different energy regions. As be- 
fore, states close to the band edge localize at small sizes 
(< 200) but the body of wavefunctions shows a much 
slower scaling behavior (compare the scales of y-axis in 
Figs. ([[5]) and ([l6|)). Since the increase of the variance 
in a given energy range is roughly the same for the four 
increases of the cluster side, a logarithmic dependence 
of the variance with the system size could be inferred 
(Notice that the number of sites is increased by a factor 
of four in each step). If this were the case, that would 
imply that only huge clusters would show exponential 
localization of wavefunctions, and consequently, Poisson 
statistics of level spacings. 
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FIG. 14. Scaling behavior of nearest neighbor spacing vari- 
ance as a function of cluster size. Squares give the evolution 
of the variance calculated for the whole band of square clus- 
ters of side L for the RHS model, circles represent the same 
magnitude for the RMF model and, filled diamonds give the 
scaling behavior for a 2 x L cluster of RMF model (in this 
case the size is not the length of the square side as before but 
\/2L) 

In conclusion, although the scaling law at the stud- 
ied sizes seems to be logarithmic for the main part of 
the band, numerical results show without doubt that 
the whole spectrum is scaling towards larger values of 
the variance, and therefore, towards localization. In any 
case, level statistics of both RHS and RMF models is not 
scaling towards the corresponding Wigner-Dyson statis- 
tics. The only alternative to a complete scaling towards 
a random distribution of levels, i.e., towards exponen- 
tial localization of the whole spectrum, is the existence 
of some new non-universal statistics for this Hamiltonian 
models. This possibility will be further discussed in the 
last Section. 
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FIG. 15. Variance of the nearest neighbor spacing distri- 
bution of the RHS model as a function of the energy plotted 
for increasing sizes of the cluster: the side of the square clus- 
ter equals 12, 24, 50, 100, and 200 for up triangles, circles, 
diamonds, squares, and down triangles, respectively. The up- 
per half of the band shows identical behavior. Error bars are 
given for the smallest and the largest cluster sizes. 
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FIG. 16. Variance of the nearest neighbor spacing distribu- 
tion of the RMF model as a function of the energy plotted for 
increasing sizes of the cluster: the side of the square cluster 
equals 12, 25, 50, 100, and 200 for up triangles, circles, dia- 
monds, squares, and down triangles, respectively. Error bars 
are given for the smallest and the largest cluster sizes. 



3. Number variance statistics 

Number variance statistics measures fluctuations in the 
number of levels that appear in an interval of fixed length. 
In this case, the statistical variable is n (number of states 
in an interval of renormalized length equal to s) and its 
variance ((n— < n >) 2 ) is calculated. Wigner-Dyson en- 
sembles produce spectra of remarkable rigidity which in 
turn means small values of the number variance |27| . On 
the other hand, a random sequence of levels gives rise 
to a number variance equal to the length of the interval. 
The advantage of using this statistics to characterize a 
level sequence is that long range correlations are explic- 
itly shown. Number variance statistic is also closely re- 
lated with the variance of the distance of one level to its 
n-th neighbor pgj . 

In accordance with the analysis done in previous sub- 
sections, number variance shows an increasing separa- 
tion of RHS (RMF) statistics from Wigner-Dyson GOE 
(GUE) values as the size of the system increases. Since it 
is really not plausible that this tendency could change be- 
yond some critical size, we can conclude that the overall 
spectra flows towards Poisson statistics, i.e., towards lo- 
calization for both RHS and RMF models. Although not 
shown in the paper, when the band is analysed grouping 
levels within different energy regions, the result of previ- 
ous subsubsection is recovered: scaling towards random 
statistics of states close to the band edges is quite rapid 
whereas scaling is much slower near band center. Once 
again, the quasi mobility edge defined by results given in 
Figs.(|]), (||) and (||) signalizes the change of behavior. 



IV. DISCUSSION 

The discussion of the results presented in Section III 
can be done according to two different scenarios. Within 
the first one, the existence of only two fixed points for the 
flow of eigenstates of a matrix belonging to one of the uni- 
versality classes is assumed. In our case, since matrices 
describing RMF model (Eq.([l])) are Hermitian (the sys- 
tem is not invariant under time reversal), the two fixed 
points are: (i) the one described by the Gaussian Uni- 
tary Ensemble (i.e., Wigner-Dyson statistics of levels se- 
quence, Porter-Thomas statistics for wavef unctions, etc.) 
and (ii), the one corresponding to localization (Poisson 
statistics for the random sequence of levels, increment 
of the number of sites having a vanishing wavefuntion 
weight, etc.). If this assumption is true, all numerical 
results prove overall scaling towards localization. More- 
over, the same trend has been proved for states belong- 
ing to different energy regions. Nevertheless, in spite of 
the clear and extense numerical evidence supporting this 
conclusion, one word of caution is necessary. Since my 
analysis is exclusively based on the scaling of averaged 
magnitudes, there is always the possibility that a count- 
able number of states of vanishing weight in averages do 
not follow overall scaling and remain extended. This con- 
clusion is in agreement with the widely accepted point of 
view that states are localized in infinite two-dimensional 
systems no matter how small the degree of disorder is [[[] . 
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SPACING 

FIG. 17. Number variance as a function of spacing for both 
RHS model (upper panel) and RMF (lower panel). Results for 
12 x 12, 24 x 24, 50 x 50, 100 x 100 and, 200 x 200 square clus- 
ters are given as continuous lines that monotonically go from 
the corresponding Gaussian ensemble result (lower dashed 
line) towards the straigt dashed line that gives the poissonian 
statistics corresponding to a random level sequence. 

On the other hand, a second scenario is also possible. 
States close to the band edges are certainly exponentially 
localized but it could happen that states within the band 
were extended, at least in a loose form: states whose spa- 
tial extension increases as P ~ L a , with a < 2, cover an 
infinite number of sites for L — > oo although the percent- 
age of visited sites relative to the total number of sites 
vanishes: P/L 2 — > 0. In this case, the observed flow to- 
wards random (localized) statistics of all measured statis- 
tical magnitudes could eventually stop at certain "new" 
fixed point. Under this assumption, statistics would start 
at Wigner-Dyson fixed point because small matrices are 
almost filled and separate from this point for larger ma- 
trices. All numerically observed effects would be assign 
to finite size effects. 

One can give further plausibility to the second option 
using results obtained for an ensemble of random real 
band matrices. It has been shown both numerically |22] ] 
and mapping the problem into a non-linear supersym- 



metric a model []29| that the average number of sites vis- 
ited by a wavefunction is ~ b 2 being b the band width. 
This gives a typical localization length of the order of b. 
In our case, we are indeed dealing with band matrices 
of a particular kind: Hamiltonian (El) on a L X L clus- 
ter of the square lattice gives band matrices of width L 
with only four nonzero matrix elements within the band. 
As a consequence, a localization length £ ~ L could be 
expected for RHS model and also, presumably, for RMF 
model. If this were the case, I think that such a large 
localization length (of the order of the sample size) can 
only be understood as signalising some kind of " geomet- 
rical" or fractal localization as the one described in the 
previous paragraph. 

Incidentally, this second scenario would be more inter- 
esting from the point of view of experimental mesoscopic 
physics than the canonical one. Instead of measuring 
universal statistics like Wigner surmise corresponding to 
the GUE, a new although non-universal statistics could 
be obtained (Compare, for example, Wigner surmise in 
the bottom panel of Fig. (|l2| ) with the distribution corre- 
sponding to 40000 sites described by the RMF model). 

In summary, there are unfortunately two different but 
internally consistent explanations of the whole set of nu- 
merical results. Either the number of fixed points is 
limited to the well known universality classes of non- 
interacting systems and then random magnetic fluxes im- 
ply localization at all energies or a new fixed point exists 
for the particular type of random band matrices ensemble 
considered in this paper and then the major part of eigen- 
states would exhibit fractal behavior instead of standard 
exponential localization. 
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